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ABSTRACT 

Context. Recent three-dimensional magnetohydrodynamical simulations have identified a disk wind by which gas materials are lost 
from the surface of a protoplanetary disk, which can significantly alter the evolution of the inner disk and the formation of terrestrial 
planets. A simultaneous description of the realistic evolution of the gaseous and solid components in a disk may provide a clue for 
solving the problem of the mass concentration of the terrestrial planets in the solar system. 

Aims. We simulate the formation of terrestrial planets from planetary embryos in a disk that evolves via magnetorotational instability 
and a disk wind. The aim is to examine the effects of a disk wind on the orbital evolution and final configuration of planetary systems. 
Methods. We perform A-body simulations of sixty 0.1 Earth-mass embryos in an evolving disk. The evolution of the gas surface 
density of the disk is tracked by solving a one-dimensional diffusion equation with a sink term that accounts for the disk wind. 
Results. We find that even in the case of a weak disk wind, the radial slope of the gas surface density of the inner disk becomes 
shallower, which slows or halts the type I migration of embryos. If the effect of the disk wind is strong, the disk profile is significantly 
altered (e.g., positive surface density gradient, inside-out evacuation), leading to outward migration of embryos inside ~ lAU. 
Conclusions. Disk winds play an essential role in terrestrial planet formation inside a few AU by changing the disk profile. In addition, 
embryos can undergo convergent migration to ~ lAU in certainly probable conditions. In such a case, the characteristic features of 
the solar system’s terrestrial planets (e.g., mass concentration around 1 AU, late giant impact) may be reproduced. 

Key words. Planets and satellites: formation - Protoplanetary disks - Planet-disk interactions 


1. Introduction 


Planets with masses larger than approximately O.lMe signifi¬ 
cantly excite density waves in a protoplanetary disk, which in 
turn exerts torque on the pla nets, and migrate toward the central 
star ( type I migration) (e.g.. lGoldreich & Tremain3ll98(]l : IWardl 
1986h . It is known that the inward migration timescale for an 
Earth-mass pl anet at 1 AU is ~ 10^ years in a locally isother¬ 
mal disk (e.g.. lTanaka et al.ll2002l) . which is shorter than the disk 
lifetime and thus can be a problem for terrestrial planet forma¬ 
tion in the solar system. Recent studies have shown that in an 
adiabatic disk, fully unsaturated positive corotation torques can 
compensate for negative Lindblad tor ques, which can change th e 
direction and rate of migration (e.g.. lPaardekooper et al.ll^lOl) . 
If this effect is significant, a convergence zone where migration 
is convergent can be cre ated in a disk and gas giant cores can 
form at a few tens of AU (iHorn et al.ll2012l) . 


iPaardekooner et aP (1201 ill revised the torque formula by 
taking into account the effects of viscous and thermal diffu¬ 
sion on the corotation torque, which indicates that only plan¬ 
ets with a limited range of masses can experience the non¬ 
linear corotation torque due to saturation effects. In fact, plan¬ 
ets with masses smaller than a few Earth masses are not affected 
by the non-linear corotation torque (the horseshoe torque) and 
they migrate inward under the influence of the Lindblad torque 
(e.g.,lKretke & Linll2012l : iHellarv & Nelsonll2012l : ICossou et al.l 

l2014l) . which means that the solar system’s terrestrial planets still 


have the problem of inward migration. However, if the slope of 
the gas surface density becomes large enough, the positive linear 
corotation torque can reverse the migration. 

During the evolution of protoplanetary disks, vertical winds 
as well as radial accretion play a significant role in the dis¬ 
persal of the gas component. While photoevaporating winds by 
UV a nd X-rays from a central star have been wide ly discussed 
(e.g., lAlexander et al.l 120061 : lErcolano et al.l 1200^ . the impor- 
tanc e of spontaneously-driven disk winds has also been pointed 
out (ISuzuki & Inutsul3120091: ISuzuki et al.ll201(il) : using mag¬ 
netohydrodynamic (MHD) simulations they showed that MHD 
turbulence triggered by magnetorotational instability (MRI) in 
disks inevitably drives vertical outflows as well as the radial 
transport of angular momentum and consequent accretion. The 
detailed properties of such di sk winds driven by turbulent Poynt- 
ing flux have been stud ied (iBai & Stone 1 120131 : lEromang et al.l 
l2013l : lLesur et al.ll2013l) . An intriguing characteristic of the disk 
wind is that disk dispersal proceeds gradually in an inside-out 
manner because the dispersal timescale is approximately propor¬ 
tional to rotation time oc where r is the radial distance from 
the star. As a result, the surface density gradient at r < lAU 
becomes much shallower than the typical values of -0.5 to -1.5, 
and could even be positive, which will drastically affect the mi¬ 
gration of protoplanets. 

As we are interested in seeing how planet formation proceeds 
in such a disk, we investigate the formation of terrestrial planets 
using A-body simulations in a disk where the effect of a disk 
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wind is considered. We calculate the disk evolution using a one- 
di mensional disk with viscous diffusion and a disk wind based 
on ISuzuki et al.l (l 2010 l) . In this work, we focus on the late phase 
of planet formation when the effect of a disk wind is clearly vis¬ 
ible, and calculate the accretion and orbital evolution of planets 
from planetary embryos with masses of 0.1 Me for ~ 100 Myr. 
The A^-body code includes the effect of type I migration that de¬ 
pends on the slope of gas surface density and the mass of the 
planet s using the prescription described by iPaardekooper et akl 
(1201 ih . We also include the influence of eccentricity on the coro¬ 
tation torque. 

Our primary aim is to clarify the effect of a disk wind on 
terrestrial planet formation. In addition, we also discuss the 
origin of the ter restrial planets in the solar system. A number 
of authors (e. g .. iKominam i & Idal 2002 : Nagasawa_et^ 2005; 
O gihara et al.l 12007 : Raymond et al.l l2009t iMorishima et al.1 
20101) have studied the origin of these planets using A^-body sim¬ 


ulations; however, no single simulation has yet convincingly re¬ 
produced the observed constraints simultaneously (e.g., eccen¬ 
tricities, spatial concentration around 1 AU, late moon-forming 
giant impact). We will propose a possible model for the origin of 
the inner solar system. 

The plan of the paper is as follows. In Sect.l^we describe our 
model and the initial conditions of the simulations; in Sect.[3]we 
present the results of A-body simulations; in Sect. 0 ] we discuss 
the origin of the terrestrial planets in the solar system; in Sect. |5] 
we give our conclusions. 


We note that although type I migration is also dependent on the 
temperature profile, we fix the temperature distribution to clarify 
the effect of the disk wind via surface density distribution. Then 
the disk scale height h is given by 


hir = 0.047 


{-ri-f 

\1AU/ \Lo/ 



( 4 ) 


where L, and M, are the luminosity and mass of the host star, 
respectively. 


2.2. Type I migration and eccentricity damping 

The orbits of planetary embryos with masses Mi, M 2 ,... and po¬ 
sition vectors ri,r 2 ,... relative to the host star are calculated by 
numerically integrating the equation of motion. 
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where Fdamp is a specific force for eccentricity and inclination 
damping and f mi„ i s a specific force for type I migration (see 
iQgihara et alJ (12014l) for each formula). The timescale for damp¬ 
ing of the eccentricity, 4 , is given by 


2. Model description 

2.1. Diskmodei 


We calculate the time evolution of the gas component of a pro¬ 
toplanetary disk by taking into account the radial mass flows by 
the transport of angular momentum through MHD turbulence 
and the mass loss due to the disk wind. Although there are un¬ 
certainties in t he mass flux of the d isk wind, we follow the simple 
framework of ISuzuki et akl (l2010l) . The evolution of gas surface 
density Sg can be expressed in the form of a diffusion equation. 


dt 


3_5 

r dr 


4/2 




■Cw^a 

V 2 ^ 


( 1 ) 


where Q is the Keplerian frequency. The standard a viscosity 
prescription, v = acji, is used, where Cs indicates the sound 
velocity and h is the disk scale height. The second term on the 
right-hand side of Eq. ([1]) is the disk wind flux pv^, which can 
be expressed as CwPoCs using t he mid-plane densi ty po, Cs, and a 
non-dimensional constant C* (ISuzuki et al.ll201(ih . We note that 
the coefficient of the first term on the rig ht-hand side of Eq. ( [T|) is 
not exactly the same as that of Eq. (9) of ISuzuki et al.l(l2010l) ow¬ 
ing to a different definition of v, which gives a slightly different 
evolution of the disk surface density. 

Eor the initial distribution of the gas disk, we assume 

S'm 2 (2) 

where /g is a dimensionless parameter and rcut = 50AU is used 
as a cutoff radius. In our A-body code, we use the evolution of 
Sg to calculate the effect of the gas disk on planet evolution. 

The tempera ture distributio n is assumed to be that of an op¬ 
tically thin disk (iHavashill 19811) . 

r.280(—) K. (3) 
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where Sg = 2400/g(r/lAU)''^^^ g cm“^ is used in the second 
equality. Here the relative motion between gas and planets is as¬ 
sumed to be subsonic (ev^ < Cs). Eor planets with high eccentric¬ 
ities and inclinations, we include a correction fac tor according to 
Eqs. (11) and (12) of ICresswell & Nelsonl(l2008h . 

The migration timescale, ta, is given by 



where is a coefficient that determ ines the direction and spee d 
of type I migration. According to IPaardekooper et akl (1201 Ol) . 
the type I migration torque depends on the Lindblad torque, the 
barotropic part of the horseshoe drag (or the non-linear corota¬ 
tion torque), the entropy-related part of the horseshoe drag, the 
barotropic part of the linear corotation torq ue, and the entropy- 
related part of the linear corotation torque. IPaardekooper et aH 
(I 2 OI ll) derived the total type I migration torque including both 
saturation and the cutoff at high viscosity. Thus we write the mi¬ 
gration coefficient in the form 


P - Pr + Pc.hmo + Pc.eaU (10) 

where/^ l, Pc.bam, andjSc.em are related to the Lindblad torque, the 
barotropic part of the corotation torque, and the entropy-related 
part of the corotation torque, respectively. Each formula is given 
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where -p and -q denote the local surface density gradient 
{p{r) - -dln'Eg/dlnr) and the local temperature gradient 
(q{r) = -d\nT/d\n r), y is the adiabatic index, and is the 
local entropy gradient (^ - q — ij - l)p)- In this work, y - 1.4 
and q = 0.5. The function F{P) is a decreasing function with 
the value of [0,1] that is related to saturation, and the functions 
G{P) and K{P) a re increasing functions wi th the value of [0,1] 
related to cutoff. iPaardekooner et akl (1201 ih introduced viscous 
and thermal parameters P^ and P^, which are expressed by the 
ratio between the viscous/thermal diffusion timescales Ty/r^- and 
the horseshoe libration timescale Tub, 



(14) 

(15) 


where the dimensionless half-width of the horseshoe region is 


Ws _ 1.1 / M r 

r y^P V Mt h’ 


(16) 


Table 1. List of parameters of each model. The variable a indicates the 
strength of turbulent viscosity while C* is a scaling factor for the disk 
wind. 


Model 

a 

Cw 

Comment 

1 

8 X 10-" 

0 

without disk wind 

2 

8 X 10-3 

2 X 10-3 

weak magnetic field 

3 

Eq. dH 

Eq. (EOll 

strong magnetic field 

4 

2 X 10-3 

5 X 10-'3 

described in Sec. 12.31 


(iFendvke & Nels^l2014il . 


^C,baro(^) 


ySc.bai-o exp 


ff 
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(17) 

(18) 


where Cf - Q.Shfr + 0.01. We note that the Lindblad torque can 
also be reduced when planets have high eccentricities and incli¬ 
nations (e.g., Papaloizou & Larwood 2000). In some runs of our 
simulations, we also add a reduction facto r to Eg . (HB, which is 
given by Eq. (13) of ICresswell & NelsonI (l2008h . However, we 
found that even if this factor is ignored, the results do not change 
very much. 

There are several works that take into account t he ef- 
fect of magne t ic fiel d on type I migrat ion ( Teryuem] 


Fromanget^ 


20051: iMuto et al.l l26o8t lUribe et al 


2003 


Bans et al. 12015 ). Among these works. iBans et al.l (12015 


2015 


con¬ 
sider the vertical transport of angular momentum by disk winds. 
In our paper, we assume the net vertical magnetic field is very 
weak and the disk wind is driven by the MRI triggered turbulent 
magnetic field, which is in contrast to the disk wind by rela - 
tively strong net magnetic field assumed in iBans et aT1 (l2015h . 
Therefore, the type I migration formula based on pure hydrody- 
na mics (which corresponds to the limit of weak magnetic field) 
by IPaardekooner et al.l (1201 Ih is probably still reasonable in our 
setting. 


2.3. Set of parameters and orbital integration 


where Ws is the half-width of the horseshoe region. For prescrip¬ 
tions of F(P), G (P), and K(P), we refer the readers to Eqs. (23), 
(30), and (31) in Paardekooner et al.l (1201 ih . 

For the thermal diffusivity, we assume = v in this work 
for simplicity. The thermal diffusivity determines the saturation 
of the entropy-related part of the corotation torque. The thermal 
diffusivity is considered to be large when the disk surface density 
decreases and the disk is optically thin, and thus radiative cool¬ 
ing is efficient. In this case, the entropy-related corotation torque 
diminishes, and one might speculate that the orbital evolution of 
planets would change. However, we do not expect that this can 
change our results in this paper. In the terrestrial pla net forming 
regio n, the disk opacity is about 1 - 10 cm^ g ' (e.g.. lBell & LinI 
Il994l) . This means that the disk is optically thin only when the 
disk is significantly dispersed (2g ~ 0.1-1 gcm“^) and type I 
migration is no longer effective. Therefore, our simplified model 
(X - v) may be justified. In addition, we observe that our qual¬ 
itative results are not affected in test simulations in which the 
entropy-related torque is switched off from the beginning of the 
simulations. 


Recent studies (e.g.. iBitsch & Kl^l2010t) suggest that the 
corotation torque decreases as the planet eccentricity increases; 
therefore, we also consider this effect using the formulae 


In our series of simulations, we mainly vary the following pa¬ 
rameters: turbulent viscosity a, disk wind efficiency Cw Al¬ 
though the efficiency of the disk wind should be related to the 
turbulent viscosity, we vary a and Cw independently to explore 
parameter space. The list of parameters for each model is sum¬ 
marized in Table [T] 

In model 1, the effect of the disk wind is not considered to 
clarify its effect on planet formation. The strength of the net ver¬ 
tical magnetic fie ld is an important control parameter that deter - 
mines a and Cw (ISuzuki et al.ll201(il: lOkuzumi & Hirosell201 ll) . 
If the net vertical magnetic field is sufficiently weak, a and Cw 
stay more or less constant during the disk evolution. Model 2 
corresponds to this weak regime; we adopt a - 8 x 10~^ and 
Cw = 2x 10~^ according to the MHD simulations of ISuzuki et al.l 
(I 2 OIOI) . If the net vertical magnetic field is relatively strong, 
a and Cw increase with time as gas materials are dispersed 
(ISuzuki et al.ll201()l) . In model 3 we consider this strong regime, 
and we use 


q: = 8 X 10 ^ X max |^1, 
Cw = 2 X 10^^ X maxll. 


0 .01Fg,i„.(r)' 

Ur) , 

0.01Fg,ini(r) 

W) 


(19) 

( 20 ) 
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In addition, in model 4 we use a = 2 x 10“^ and = 5 x 10^^. 

The parameter set for model 4 may require some explanation. 
The slope of gas surface density basically depends on a/Cw. In 
other words, if the ratio of a/Cw is the same, the distribution 
of gas surface density takes the same form while the evolution 
timescale depends on a and Cw. In model 4 we consider the case 
of a more efficient disk wind (a/Cw = 40) than that of model 2 
(a/Cw = 400). The value of a is also changed because we wish 
to discuss the origin of the solar system’s terrestrial planets in 
this paper. As we show in Sects. 0 and 01 outward migration 
of 0.1 Earth-mass embryos is favorable for reproducing the ob¬ 
served properties of the solar system and a ~ 10 “^ is required 
to avoid the saturation of the corotation torque for these bodies. 
Therefore, a -2x 10“^ and Cw = Q'/40 = 5 x 10“^ are adopted. 
While the evolution of the surface density of model 4 is slower 
than models 2 and 3 because both a and Cw are smaller, the trend 
of the radial dependence of the surface density of this model 
represents a case more or less between model 2 and model 3. 
The initial amount of gas and the gas dispersal time are also 
changed in model 4. The scaling factor for the initial gas den¬ 
sity /g is increased by a factor of three. In the other models, 
/g is set to unity. In addition, the disk would rapidly disperse 
(~ O.lMyr) after typical disk lifetimes of a few Myr, and addi¬ 
tional effects ar e required to achieve t his (e.g., photoevaporat¬ 
ing winds) (e.g.. [Alexander et al.ll2014)l . In model 4, in order to 
mimic this disk evolution, a and Cw are increased by a factor of 
100 after t = 3Myr. 

For the initial distribution of embryos, we put 60 protoplan¬ 
ets with mass M - O.lMe between a = 0.1 - 3 AU. Their or¬ 
bital separation. A, measured in mutual Hill radii is set to 10. We 
perform three runs of simulations for each model with different 
initial locations of protoplanets. 

For numerical integration, we use a fourth-order Hermite 
scheme (iMakino & AarsetBll992h . When the physical radii of 
two bodies overlap, they are assumed to merge into one body, 
conserving total mass and momentum assuming perfect accre¬ 
tion. The physical radius of an embryo R is determined by 
R - ( 3 M/ 4 ; 7 rp)*^^, where we assume the internal density of 
p - 3g cm~^. 

3. Results 

3.1. Disk evolution 

We first show the disk evolution and its contribution to type I 
migration. The solid lines in Fig. [T] show the evolution of the 
disk surface density for models 2-4 while the dotted lines in each 
panel represent the surface density for model 1. Figure |2] shows 
the migration timescale (C) for model 2 and model 4, where the 
migration rate and direction are indicated by the color (see color 
bar). When the migration timescale is negative, the direction of 
migration is inward. As we saw in Sect. |2l type I migration de¬ 
pends on the eccentricity of planets. We assume e = 0.01 in these 
plots. 

Figure [n a) clearly shows that the disk evolution of model 2 
is altered by the disk wind in the following two aspects. First, 
the slope of the disk becomes shallow in the inner region. Sec¬ 
ond, the disk is depleted faster than the case without the disk 
wind. These two factors significantly reduce the inward migra¬ 
tion of embryos. From the migration map of model 2 in Fig.|3a), 
we find that the inward migration timescale for sub-Earth mass 
planets is over 1 Myr at f = O.lMyr, which means they no longer 
undergo migration. We also observe no migration at f = IMyr 
and 3 Myr. 



Fig. 1. Evolution of gas surface density profile for (a) model 2, (b) 
model 3, and (c) model 4. Solid lines in each panel indicate the gas 
surface density at r = 0,O.lMyr, IMyr, and 3Myr, respectively. Dotted 
lines show results for model 1. In panels (a) and (b), the solid and dotted 
lines overlap for r = 0 because of the initial conditions. 

For the disk evolution of model 3 that has larger C* and a 
shown in Fig.lHb), an inner cavity is created inside a certain ra¬ 
dius at which 99% of the initial gas is dispersed according to Eqs. 
(fT^ and (l20l) . As time passes, the inner cavity grows gradually 
and the disk inner edge moves outward. The inside-out evacua¬ 
tion of the disk can significantly affect the orbital evolution of 
planets, which will be shown in Sect. 13.21 

We also examine the disk evolution for model 4, in which 
the radial profile of the surface density in the inner region is 
approximately between model 2 and model 3. In Figs. [T]and| 2 ] 
for model 2 , we see the disk is depleted in the early stage of 
planet formation; therefore, a larger amount of initial gas disk 
is assumed (fg - 3). The evolution of the surface density in 
Fig- Ec) shows that the slope becomes positive inside 1 AU. 
In order to get a positive surface density gradient inside 1 AU, 
we find that a/C* ^ 100 is required. The disk lifetime is longer 
than model 2 , which is more consistent with the observation- 
ally inferred disk lifetime. The migration map in Fig. |2^b) is 
more complicated. Because the effect of the disk wind is stronger 
in the inner region, the orbital region where planets can move 
outward is limited to r < lAU. The saturation of the coro¬ 
tation torque depends on the mass and the disk viscosity. If 
a = 2 X 10“^ is assumed, the saturation effect is weakest for 
bodies with M ^ 0.1 M^. Therefore, only sub-Earth-mass plan¬ 
ets can undergo outward migration. It is also worth noting that 
planetary orbits can evolve until a few Myr because the disk is 
almost completely depleted after 10 Myr. 

Before proceeding to the next subsection, where the results 
of A-body simulations are shown, here we briefly discuss how 
the inner disk profile is determined by viscous diffusion and the 
disk wind. We note that the disk wind brings about the variety in 
the radial dependence of the surface density in the inner region 
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Fig. 2. Time evolution of migration rate and direction for bodies with e = 0.01 for (a) model 2 and (b) model 4. The color indicates the migration 
timescale. The contours show migration timescales of 0.1 Myr, 1 Myr, and 3 Myr. 


ranging from a clear hole (model 3) to mild inside-out evacuation 
(model 2). We suppose that the inner hole seen in model 3 is an 
extreme case. In a realistic situation low-density gas should re¬ 
main in the inner region and a low level of accretion would con¬ 
tinue, rather than the formation of a clear inner hole, even at later 
times as seen in model 4, because the origin of the energy to drive 
the d isk wind is the grav itational energy liberated by gas accre¬ 
tion (ISuzuki et al.ll20ldl) : the infall of a certain amount of gas is 
required to launch vertical outflows. ISuzuki et al.l(l2010h demon¬ 
strated the time evolution of a global disk by a simple model tak¬ 
ing into account such energy-limited disk winds, which shows 
that the radial profile of the surface density is slightly modified 
at < 1 AU. 

3.2. Orbital evolution 

We show the results of A-body simulations in disks that evolved 
via viscous diffusion and a disk wind. Figure[2shows the orbital 
evolution of planets. The size of each filled circle represents the 
radii of bodies. The color of lines indicates the eccentricity of 
the planets. Each panel shows the results for model 1, model 2, 
model 3, and model 4. We performed three runs for each model, 
and typical runs are shown in the plots. Figure |4] shows the final 
orbital configurations in the semimajor axis-mass plane for all 
runs of model 2, model 3, and model 4 (the results for model 1 
are not shown), where planets that formed through the same run 
are connected with the same line. 

In Fig. [3a), the results for model 1 (in which the disk wind is 
not considered) is shown for comparison. Planetary embryos mi¬ 
grate inward and 12 bodies are lost by f = 20Myr. The computa¬ 
tional cost is huge for this run because the bodies do not collide 
with each other and the orbital period of the innermost planet 
is short; therefore, the computation was stopped at f = 20Myr, 





Radial Distance (AU) 

Fig. 4. Final orbital configurations for (a) model 2, (b) model 3, and (c) 
model 4. Only for model 4a is the configuration at f = 200Myr shown. 
Results of three different runs are plotted in each panel. 

which took over four months of CPU time. Orbit crossings are 
not observed before t - 20Myr; however, the system is expected 
to be destabilized at ~ 100 Myr acco rding to previous studi es on 
the stability of multi-planet systems (IChambers et alJ[T996h . 

In the results for model 2, we find that planetary embryos 
do not migrate inward any more as expected from Figs. [1] and 
13 In this model, as the disk dissipates, the eccentricity of the 
embryos is excited and they exhibit orbit crossings, leading to 
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Fig. 3. Orbital evolution of planetary embryos for (a) model 1, (b) model 2, (c) model 3, and (d) model 4. The filled circles connected with lines 
represent each body. The color of lines indicates the eccentricity (see color bar). The size of circles is proportional to the radius of the body. 


giant impacts between the bodies. At f = lOOMyr, the number of 
planets has decreased to 13 and their eccentricities are between 
about 0.01 and 0.1. The final orbital configurations are almost the 
same among the three different runs (see Fig. 01 a)). The averaged 
mass of the largest body in each run is ^ O.OMe. The typical 
orbital separation between the planets is ^ 10 - 40rH. 

For model 3, the orbital evolution is significantly different 
from that without a disk wind. As seen in Fig. 01b), the disk’s 
inner edge moves outward. The migrating inner edge sweeps up 
embryos and almost all bodies migrate outward. As a result, no 
planets are left inside 0.4 AU and a large planet forms at 5 AU. 
The eccentricities of the planets are ~ 0.04-0.3 (some bodies are 
still exhibiting orbit crossings). In the other two runs, the results 
are almost the same (see Fig. 01b)). The mass of the largest body 
is ^ 3.8Mo, while the second largest planet is located at around 
lAU with a mass of ^ 0.5M®. 

We observe an interesting orbital evolution in model 3 in 
which most protoplanets are swept up by the disk edge before 
t ^ 0.3Myr but afterward several bodies drift through the edge. 
We briefly discuss how bodies follow the migration of the disk 
edge or not. First, just as there are for torques on planets, there 


are two different positive torques, namely the type I migration 
torque and the edge torque. As seen in Sect.|2] the former can be¬ 
come positive owing to the corotation torque, which depends on 
the slope of the gas surface den sity. The latter is a torque recently 
found bv lQgihara et al.l (l2010l) . When a body with non-zero ec¬ 
centricity straddles a sharp disk inner edge, the body gains a 
net positive torque from the disk depending on the sharpness 
of the edge and on the eccentricity. Whether planets can follow 
the movement of the disk edge depends on the migration speeds 
of the planets and the edge. If the outward migration rate of the 
planets due to the positive torque (which is the sum of the type 
I migration torque and the edge torque) is smaller than the mi¬ 
gration rate of the edge, the planets go through the disk edge 
and stay in place. Through a series of test calculations, we find 
that the positive torque exerted on the planets is usually large 
for this model and 0.1 Earth-mass planets can follow the disk 
edge. Therefore, there is another reason why some planets drift 
through the disk, which is that small embryos are scattered in¬ 
side the disk edge by larger bodies. Larger bodies grow at the 
disk edge by sweeping up smaller embryos. At f ^ IMyr, the 
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planets at the edge already have masses of ^ 3.5M® and scatter 
embryos inside the edge. 

The results for model 4 are also interesting. Before t ^ 0.02 
Myr, all bodies undergo inward migration, but afterward ~ 0.1 
Earth-mass embryos in close-in orbits move outward while bod¬ 
ies in distant regions migrate inward, leading to a convergence 
of planets at around 1 AU. For a - 2 x 10“^, the positive 
corotation torque can be saturated for bodies larger than ^ 0.3 
Earth masses; as a result, several planets that undergo collisional 
growth migrate inward and are lost to the central star. As the 
disk disperses after t - 3Myr, all the bodies cease their migra¬ 
tion. At f = lOOMyr, a system with six planets with a maximum 
mass of Q.lM(s and e ^ 0.01 had formed. The orbital separations 
between the planets are relatively small (down to < lOrn), and 
giant impacts can occur after t - lOOMyr. Thus we continued 
simulation until t - 200Myr for this run. The system underwent 
a giant impact at f = 145Myr, which is discussed in Sect. H] In 
the other two runs, although the mass distributions are slightly 
different from each other, the results are qualitatively the same 
(Fig. Etc)). 


Table 2. Statistical measures for final planetary systems. Sc is a measure 
of the radial mass concentration. Sd measures how a system deviates 
from circular and coplanar orbits. The last column shows the timing of 
the last giant impact, although last impact events do not occur within 
the simulation time for model 1 and model 4 (except for 4a). 


Run 

5c/89.9 

5d/0.0018 

^impjast (Myr) 

la 

0.042 

0.0016 

> 20 

lb 

0.050 

0.0016 

> 20 

Ic 

0.050 

0.0013 

> 20 

2a 

0.048 

1.7 

61 

2b 

0.052 

1.5 

30 

2c 

0.047 

9.1 

80 

3a 

0.086 

3.7 

99 

3b 

0.12 

1.5 

57 

3c 

0.098 

0.87 

72 

4a 

0.80 

0.046 

145 

4b 

0.81 

0.0045 

> 100 

4c 

0.35 

0.020 

> 100 


4. Origin of the soiar system’s terrestriai pianets 


Here we discuss our results in terms of the origin of terrestrial 
planets in the solar system. The important characteristics are 
the radial mass concentration and small eccentricities of planets, 
which we characterize using two statistica l measures. Sr. an d Sa. 
The radial mass concentration statistic is (IChambersll200 Ih 


Sc - max 






( 21 ) 


where j - 1,2, The angular momentum deficit, wh ich rep¬ 
resen ts the deviation from circular and coplanar orbits, is (iLas 
l1997h 


2 nij 




■ cosiij) 


2 


( 22 ) 


where ij is the orbital inclination of the body j. 

The radial mass concentration of the solar system’s terrestrial 
planets is = 89.9, which indicates high radial mass concentra¬ 
tion because the masses of Mercury and Mars are small and the 
orbital distance between Venus and the Earth is also relatively 
small. Previous studies have shown that it is difficult to repro¬ 
duce such a high mass concentration te.g.. lRavmond et ^120091: 
iMorishima et al.ll2010l) : a way to account for this statistic in N- 
body simulations is to ass ume a mass distribution initially con - 
fined to a narrow annulus (IMorishima et al.l2008UHanse nl20n9h . 
or in other w ords to assume an in i tially large 5 p. In the grand tack 
model (e.g., IWalsh et"^ 1201 it iJacobson & Morbidellil 120141) . 
Jupiter sweeps up to ~ 1.5AU and embryos are eventually con¬ 
fined to a relatively compact region. However, they assumed an 
inner cutoff at ^ 0.7AU, which means that is already rela¬ 
tively large in the initial state. 

The angular momentum deficit for the terrestrial planets in 
the solar system is 5d = 0.0018, which means that the eccen¬ 
tricities of the terrestrial planets in the solar system are relatively 
small (e ^ 0.01 - 0.1). We note that Jupiter and Saturn increase 
S d after their formation, which is ignored in this paper. There¬ 
fore, the necessary condition for our results to explain terrestrial 
planets in the solar system is that Sd < 0.0018. 

Table |2] shows the statistical measures for all runs. In the 
cases of model 1 (la, lb, Ic) and model 2 (2a, 2b, 2c), Sc is kept 


small from the beginning of the simulations because the planets 
do not undergo convergent migration. The initial value of Sc,o is 
5.4. In the results of model 3 (3a, 3b, 3c), is slightly increased 
to ^ 9. The final system is more concentrated than the initial dis¬ 
tribution, but the value is still not consistent with the value of the 
current solar system. In the results of model 4 (4a, 4b, 4c), 5c is 
as large as the value of the inner solar system {Sc is between 30 
and 70). The mass concentration at around 1 AU is clearly seen 
inFig.[3ld). 

Regarding the angular momentum deficit, Sd for model 1 
is small because planets do not exhibit orbit crossings. In con¬ 
trast, in model 2 and model 3, planets experience orbit crossings 
and the gas disk is depleted early, thus 5 d is larger than that of 
model 1. We neglect the forced eccentricity by Jupiter and Sat¬ 
urn such that additional damping may be required for model 2 
and model 3. In the results of model 4, 5d is sufficiently small 
(r^ 4 X 10-5). 

The timing of a moon-forming impact is also a key con¬ 
straint on formation models of the solar system. According to 
recent work, the time of the last giant impact on Earth could 
have been significantly late (50-150 Myr a fter the formation of 
the first solids in the solar systern) (e.g ., iTouboul et al.l 120071 
lAllegre et al.ll200^ IJacobson et al.ll201^ . The timing of the lat¬ 
est impact for model 2 is f ~ 50 Myr, and that for model 3 is 
f ^ 50 - 100 Myr, consistent with the timing of the last giant 
impact on Earth. Particularly interesting is that the timing of the 
last giant impacts is delayed after 100 Myr for model 4. In fact, 
we continue simulation for 4a until t = 200Myr and find that the 
last giant impact occurs at f = 145Myr. 

In conclusion, we find that a disk wind can increase the ra¬ 
dial mass concentration. In particular. Sc can be as large as that 
of the solar system even if we assume that the initial embryo dis¬ 
tribution is extended to 0.1 AU. In the result of 4a, the orbital 
properties including the timing of the moon-forming impact are 
similar to those of the solar system. We have only performed 
simulations in a limited number of cases; however, we can iden¬ 
tify a significantly large parameter space that may allow us to 
find a best-fit model for the observed properties of the solar sys¬ 
tem. To assess the possibility, we have to explore a broad range 
of parameter space to find a condition that can create a solar sys¬ 
tem analog. 

We also discuss one other problem in the formation of the 
solar system. Recent numerical simulations on the evolution of 
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the snow line in the solar system that solve the detailed radiative 
energy transfer have shown that the snow line moves in the disk 
and comes inside Earth’s orbit to ^ 0.5AU (e.g. JOka et al.ll20Tll : 
iDavisIImpost iGraud & Linll200'^ during a phase with a disk ac¬ 
cretion rate of ~ 10 '®Moyr If planetesimals formed during 
this stage, the terrestrial planets in the solar system should have 
formed from icy planetesimals and contain significant amounts 
of water, which is not consistent with the current water content 
of the terrestrial planets. We can provide a pathway to avoid this 
problem. Protoplanets that formed from rocky planetesimals in¬ 
side the snow line within < 0.5 AU migrate outward in disks that 
are affected by a disk wind, resulting in water-poor planets at 
around 1 AU. 

5. Conclusions 

We performed A-body simulations that include the effects of a 
disk wind and found that the orbital evolution of terrestrial plan¬ 
ets and their precursor protoplanets can significantly differ from 
the results of previous simulations. Because of the disk wind, the 
efficiency of which depends approximately on gas materi¬ 
als blow out from the surface of the disk from the inside out, 
leading to a shallower surface density profile in the inner re¬ 
gion (< lAU). Planets with masses less than a few Earth masses 
would no longer migrate inward. The orbital evolution depends 
on the parameters (a, €„), and in fact, 0.1 Earth-mass embryos 
inside 1 AU never undergo significant inward migration for all 
the parameters in our simulations. In an extreme case in which 
the vertical magnetic field is strong, we observed the inner edge 
of the disk move outward and sweep up embryos. 

In a certain parameter range, we found that planetary em¬ 
bryos move outward. When a/C* is less than ~ 100, the surface 
density slope becomes positive in the inner region, which can in¬ 
duce convergent migration to around 1 AU. In order for outward 
migration to take place, the corotation torque should not be satu¬ 
rated, which depends on diffusion. Eor 0.1 Earth-mass embryos, 
the saturation is inhibited when a is ^ 10“^. Thus, we performed 
additional simulations with a - 2 x 10“^ and C* = 5 x 10“^ 
and found that the resultant systems are radially concentrated at 
around 1 AU. This is an important result. In this case, the ob¬ 
served constraints of the inner solar system (e.g., eccentricity, 
radial mass concentration, late moon-forming impact) may be 
reproduced. This expectation requires that the turbulent viscosity 
be quite low (a ~ 10“^) in the terrestrial planet-forming region 
in the protosolar nebula. Although we focus on the properties 
of the terrestrial planets, other observed properties in the solar 
system (e.g., the dynamically excited asteroid belt) should be re¬ 
produced in future work. A combinat ion of our model with the 
grand tack model (IWalsh et al.ll201 ih . which explains asteroid 
excitation, may be a possible scenario. In addition, a disk wind 
can play other important roles. Eor example, if the snow line had 
been located inside 1 AU, the Earth would have accreted a sig¬ 
nificant amount of icy planetesimals. However, this problem can 
be overcome if the Earth formed from rocky protoplanets that 
migrated outward owing to the effect of the disk wind. 

The disk wind can also apply to exoplanet formation models. 
A large number of close-in low-mass planets have been discov¬ 
ered, and thus we can statistically compare observational results 
with theoretical models. To do so, we need to explore parameter 
space (a, Cw), which we leave for future work. If disk winds op¬ 
erate efficiently in exoplanet systems, super-Earths would tend 
to pile up at around 1 AU rather than 0.1 AU. 
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